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Abstract 

We consider sustainable management issues formulated within the framework of 
control theory. The problem is one of controlling a discrete-time dynamical system 
(e. g. population model) in the presence of state and control constraints, representing 
conflicting economic and ecological issues for instance. The viability kernel is known 
to play a basic role for the analysis of such problems and the design of viable control 
feedbacks, but its computation is not an easy task in general. We study the viability 
of nonlinear generic ecosystem models under preservation and production constraints. 
Under simple conditions on the growth rates at the boundary constraints, we provide 
an explicit description of the viability kernel. A numerical illustration is given for the 
hake-anchovy couple in the Peruvian upwelling ecosystem. 
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1 Introduction 



This paper deals with the control of discrete-time dynamical systems of the form 
x(t + 1) = f (x(t),u(t)j , t G N, with state x(t) G X and control u(t) G U, in the presence 
of state and control constraints (x(t),u(t)) G D. The subset D C X x U describes "ac- 
ceptable configur ations o f the s ystem" . Such p roblems of dynamic control under cons t raints 
refer to viability (lAubinl . Il99ll ) or invariance (IClarke. Ledayev. Stern, and Wolenskil . Il995l ) 
frameworks. From the mathematical viewpoint, most of viability and weak invariance results 
are addressed in the continuous time case. However, some mathematical works deal with 



the discrete-time case. This includes the study of numerical s chemes 



tion o 



1994 



the viability problems of the c ontin uous dynamics as in ( lAubin 



or the approxima 



1991 



Saint-Pierrel . 



Quincampoix and Saint-Pierrd . Il995l ). In the control theory literature, problems of 



constrained cont rol have also been addressed in the discrete time case (see the survey paper 



(IBlanchini 



1999j0; reachability of target set s or tubes for nonlinear discrete time dynamics 



is examined in (IBertsekas and Rhodes 



197lh 



We con sider sustainable managemen t issues which can be formu 



work as in ( 



Bene and Doven 



2000 



2003 



Bene. Doven. and Gabav 



atec. 



2001 



within such a frame- 



Eisenack. Sheffran. and Kropp 



2006; M artinet and Dovenl. 120071; iMullon. Cury and Shannon! 12004; iRapaport. Terreaux. and Doyen 



2006 



De Lara. Doyen. Guilbaud. and Rochet 



2007 



De L and Dovenl 



2008 



The time index t is an integer and the time period [t, t + 1[ may be a year, a month, 
etc. The dynamics is generally a population dynamics, with state vector x(t) being either 
the biomass of a single species, or a couple of biomasses for a predator-prey system, or 
a vector of abundances at ages for one or for several species, or abundances at different 
spatial patches, etc. The control u(t) may represent catches, harvesting mortality or har- 
vesting effort. The set P may i nclude biological, ecological and economic objectives as in 



( iBene. Doyen, and Gabayl . 1200 ll ). For instance, if the state x is a vector of abundances at 
ages and the control u is a harvesting effort, D = | B(x) > 6 b ,E(x,u) > e 9 } repre- 

sents acceptable or acceptable configurations where conservation is ensured by an biological 
indicator B(x) > $ (spawning stock biomass above a reference point, for instance) and 
economics is taken into account via minimal catches E(x, u) > e° (catches E(x, u) above a 
threshold). 

The viability kernel V(/, ©) associated with the dynamics / and the acceptable set D is 
known to play a basic role for the analysis of such problems and the design of viable control 
feedbacks. Unfortunately, its computation is not an easy task in general. In this paper, we fo- 
cus on the case where the viability kernel is given by a finite algorithm, with emphasis on pop- 
ulation dynamics ecosystem models. For nonlinear predat or-prey systems, one can find de- 



scriptions of viability kernels for Lotka-Volterra systems in (IBonneuil and Miillers 



popu lation viability analysis in three trophic-level food chains in (IBonneuil and Saint-Pierre 



1997) and 



20051 ). In these latter references, the control strategies are carried by the species ("fitness- 
maximizing" strategies), the constraint set corresponds to having high enough densities, and 
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the models are in continuous time. In our approach, the controls are the harvesting efforts, 
the constraint set includes both high abundances and high catches, and the model is in 
discrete time. 

The paper is organized as follows. Section [2] is devoted to recalls on discrete-time viability 
and its possible use for sustainable management. For generic nonlinear ecosystem models, 
we provide an explicit description of the viability kernel for preservation and production 
constraints in Section [21 together with viable controls. An illustration in ecosystem man- 
agement and numerical applications are given for the hake-anchovy couple in the Peruvian 
upwelling ecosystem in Section HI 

2 Discrete— time viability 

Let us consider a nonlinear control system described in discrete-time by the difference equa- 
tion 

x(t + 1) = f{x(t), u(t)) for all t E N, 
x(Q) = xq given, 

where the state variable x(t) belongs to the finite dimensional state space X = M nx , the 
control variable u(t) is an element of the control set U = M nu while the dynamics f maps 
X x U into X. 

A controller or a decision maker describes "acceptable configurations of the system" 
through a set D C X x U termed the acceptable set 

(x(t),u{t)) E D for all t E N , (2) 

where D includes both system states and controls constraints. 

The state constraints set V° associated with D is obtained by projecting the acceptable 
set D onto the state space X: 

V° : = Proj x (D) = {x E X | 3u E U , (x, u) E D} . (3) 

Viability is defined as the ability to choose, at each time step t E N, a control u(t) E U 
such that the system configuration remains acceptable. More precisely, the system is viable 
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if the following feasible set is not empty: 



V(/, 



Xq E X 



3 («(0),«(1),...) and (x(0),x(l),. 
satisfying (pP) and ((2J) 



(4) 



The set V(/, D) is called the viability kernel (jAubinl . Il99ll ) associated with the dynamics / 



and the acceptable set D. By definition, we have V(/, D) C V° = Proj x (D) but, in general, 
the inclusion is strict. For a decision maker or control designer, knowing the viability kernel 
has practical interest since it describes the states from which controls can be found that 
maintain the system in an acceptable configuration forever. However, computing this kernel 
is not an easy task in general. 

We now focus on the tools to achieve viability. A subset V is said to be weakly invariant 
for the dynamics / in the acceptable set D, or a viability domain of / in D, if 



VxGV, 3ueU, (x,u) E D and f(x,u) E V . 



(5) 



That is, if one starts from V, an acceptable contro l 



Moreover, according to viability theory (lAubin 



may transfer the state in V. 



19911 ). the viability kernel V(/,D) turns 



out to be the union of all viability domains, or also the largest viability domain: 



V(/, 



IjJ V, V C V°, V viability domain for / in 



(6) 



A major interest of such a property lies in the fact that any viability domain for the dynamics 
/ in the acceptable set D provides a lower approximation of the viability kernel. An upper 
approximation of the viability kernel is given by the so called viability kernel until time 
k associated with f in D: 



Y, := { 



Xq E X 



3 (u{0),u{l),...,u{k)) and (z(0), x{l), . . . , x{k)) 
satisfying (jTJ) for t — 0, . . . , k — 1 



> . 



and d25 for t — 0,...,k 



(7) 



We have 



V(/, D) C V fc+ i C V fc C V = V° for all k E N . 



It may be seen by induction that the decreasing sequence of viability kernels until time k 
satisfies 



V = V° and Y k+1 = {x E Y k \ 3u E U , (x, m)g1 and f(x, u) E W k } . 



(9) 



By OH]), such an algorithm provides approximation from above of the viability kernel as 
follows: 



V(/, 



lim IV* 



(10) 



Conditions ensuring that equality holds may be found in (ISaint-Pierrd . Il994l ). Notice that, 
when the decreasing sequence (Vfc)feeN of viability kernels up to time k is stationary, its limit 
is the viability kernel. Indeed, if = V k+ i for some k, then is a viability domain by (Q. 
Now, by dlj), V(/, D) is the largest of viability domains. As a consequence, V k = V(/, O) 
since V(/, D) C V* by ((SJ). We shall use this property in the following Sect. [31 

Once the viability kernel, or any approximation, or a viability domain is known, we have 
to consider the management or control issue, that is the problem of selecting suitable controls 
at each time step. For any viability domain V and any state x G V, the following subset 
Uv(x) of the decision set U is not empty: 



Uv(x) := {u E U | (x,u) E © and f{x,u) E V} . 



11) 



Therefore Uv(/,d)(^) stands for the largest set of viable controls associated with x E X. Then, 
the decision design consists in the choice of a viable feedback control, namely any selection 
: X — ► U which associates with each state x E V(/, D) a control u = ^>(x) satisfying 
#0) E U V (/,D)(ar). 

In the context of sustainable management, viability concepts and methods may help giv- 
ing a framework for setting decision making. First, one should delineate perpetual objectives 



approach is illustrated in ( 


De Lara. Doven. Guilbaud. and Rochet. 


2006 


), and especially in 


(De Lara. Doven. Guilbaud. and Rochet. 


2007 


) for fishery management. 
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3 Viable control of generic nonlinear ecosystem models 



For a generic ecosystem model, we provide an explicit description of the viability kernel. 
Then, we shall specify the results for predator-prey systems, in particular for discrete-time 
Lotka-Volterra models. 



3.1 Dynamics and constraints for an ecosystem management model 

For simplicity, we consider a two-dimensional state model. However, the following Proposi- 
tion [1] may be easiliy extended to n-dimensional systems as long as each species is harvested 
by a specific device: one species, one harvesting effort. The two-dimensional state vector 
(y, z) represents biomasses and the two-dimensional control (v, w) is harvesting effort of each 
species, respectively. The catches are thus vy and wz (measured in biomass)Jl 
The dynamics / in ([I]) is given by 

(yR v (y,z,v) \ 
for all (y, z, v, w) £ I , (12) 
zR z {y,z,w) J 

where R y : R 3 -> R and R z : M 3 -> R are two functions representing growth coefficients (the 
growth rates being R y — 1 and R z — 1). We do not make assumptions on the signs of R y and 
R z ; the restrictions on the domain of variation of (y, z, v, w) will result from the requirement 
that the trajectories belong to the acceptable set D, which will include signs considerations. 



The discrete-time control system ([T]) with dynamics / given by (fl2l now writes: 

y(t+l) = y(t)R y (y(t),z(t),v(t)) 
z(t + l) = z(t)R z (y(t),z(t),w(t)) . 

The implicit assumption in such functional form (fl2j) for the dynamics is that, during one 
time period, the harvesting effort v(t) of species y only affects the same species y(t + 1) and 



fact, any expression of the form c(y, v) would fit for the catches in the following Proposition[T]as soon 
as v i— > c{y, v) is strictly increasing and goes from to +oo when v goes from to +oo. The same holds for 
d(z, w) instead of wz. 
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not the other one zit + 1) (the reverse holds for w(t) and species z). Of course, after two 
periods, zit + 2) depends on y(t + 1) which depends on v(t) so that both efforts affect both 
species. Thus, the time period is assumed to be short enough for the impact on one species 
of harvesting the other species to be negligible. 

The acceptable set D is defined by minimal biomass thresholds y 9 > 0, z 9 > and minimal 
catch thresholds Y 9 > 0, Z 9 > 0: 

D = { (y, z, v, w) E R 4 | y > y\ z > z\ vy > Y\ wz>Z h }. (14) 

Such a set includes both preservation 

y(t) > y b , z(t) > z b , 

and production requirements 

v{t)y{t) > Y b , w(t)z(t) > Z b 

when thresholds are positive. 

3.2 Expression of the viability kernel 

The following Proposition [T] gives an explicit description of the viability kernel, under some 
conditions on the minimal thresholds. 

Proposition 1 Assume that the function R y : IR 3 — > R is continuously decreasing in the 
control v and satisfies lim v ^ +00 R y (y , z , v) < 0, and that R z : R 3 — ► R is continuously 
decreasing in the control variable w, and satisfies \im. w ^ +oa R z (y, z, w) < 0. If the thresholds 
in (TT41) are such that the following growth coefficients are greater than one 

R y {y\ z\ > 1 and R z {y\ z\ \)>\, (15) 
y" z" 

2 In all that follows, a mapping ip : R — > R is said to be increasing if x > x' => > fi^)- The reverse 
holds for decreasing. Thus, with this definition, a constant mapping is both increasing and decreasing. 
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the viability kernel associated with the dynamics / in (|1~2|) and the acceptable set 
is given by 

-b 



in (EE 



v(/,: 



{yb p\> 
{y, z)\y>y\ z> z\ yR y (y, z, —) > y\ zR z (y, z,—)>z b j. 



(16) 



Before giving the proof, let us comment the assumptions. That the growth coefficients 
are decreasing with respect to the harvesting effort is a natural assumption. Conditions (fl~5l) 
mean that, at the point (y b , z b ) and applying efforts uf = v b = the growth coefficients 
are greater than one, hence both populations grow; hence, it could be thought that computing 
the viability kernel is useless since everything looks fine. However, if all is fine at the point 
(y\z b ), it is not obvious that this also goes for a larger domain. Indeed, the ecosystem 
dynamics / given by (1121) has no monotonocity properties that would allow to extend a 
result valid for a point to a whole domain. What is more, if continuous-time viability results 
mostly relies upon assumptions at the frontier of the constraints set, this is no longer true 
for discrete-time viability. 

Proof. According to induction Q, we have: 



Vn 



Vi 



{(y,z) y >y ,z>z }, 



{v,z) 



(y,z) 



y > y", z > z b and, for some (v, w) > 0, 
vy > Y\wz > Z\ yR y (y, z, v) > y\ zR z (y, z, to) > z b 



y>y,z>z, yR y (y, z, — ) > y , zR z (y, z, — )>z 



y 



because v i— > R y (y, z,v) and w i— > R z (y, z,w) are decreasing, 

yb Z \> 

and thus we may select v = — , w = — . 

y z 

Denoting y' = yR y (y,z,v), z' = zR z (y, z,w), we obtain, 



V, 



{y,z) 



y > y° , z > z and, for some (v, w) > 0, 
vy > Y\ wz > Z b 

y">y\ y'Ry(y',z'X)>y\ z">z\ z> R z (y' , z> , f) > z b 



We shall now make use of the property, recalled in Sect. [21 that when the decreasing sequence 
(Vfc)fceN °f viability kernels up to time k is stationary, its limit is the viability kernel V(/, B). 
Hence, it suffices to show that Vi C V2 to obtain that V(/, B) = Vi. Let (y, z) E Vi, so that 



yb ff> 
y>y b , z > z b and yR y (y,z, — ) > y, zR z (y,z, — ) > z b . 

y z 



Since R y : R 3 — > R is continuously decreasing in the control variable, with lim,_,_> +00 R y (y, z, v) < 0, 
and since yRy(y, z,^-) > y b , there exists a v > ^- (depending on y and z) such that y' = 
yR y (y, z, is) = y" . The same holds for R z : R 3 — > R and z' = zR z (y, z, w) = z b . By (fT5j) . we deduce 
that 

y'i^S/V,-^) =y b i2 y (y b ,^X) >y b and z'it^y', ^) = z b R z (y b , z\ %) > z b . 
y y z z 

The inclusion Vi C V2 follows. □ 

As a direct consequence of the proof of Proposition fTJ for each (y, z) e V(/, D), the control 
(v, w) belongs to Uv(/,b)(2/> z ), the set of viable controls denned in (ITTT) . More explicitly, 
we have the following expression of the viable controls set, which results from the above 
observation and proof. 

Corollary 2 Suppose that the assumptions of Proposition^ are satisfied. Denoting 

Y b 

v(y, z) = max{v > — \ yR y (y, z, v) = y b } 

Z b 

w(y,z) = max{u> > — | zR z (y, z,w) = z b } , 



the set of viable controls is given by 
where y 1 = yR y (y,z,v), z' = zR z (y,z,w). 



v(y, z) > v > — , w[y, z) > w > — , 
y z 

y'Ry(y', z', ^j) > y\ z'R z (y', z\ > z b 
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3.3 Viability kernel for a nonlinear predator— prey system 

Till now, no trophic relationship has been specified between species y and z. From now 
on, we shall focus on nonlinear predator-prey systems. The two-dimensional state vector 
(y, z) represents biomasses of preys (y) and predators (z). This interpretation results from 
the monotonicity assumptions made on function R z in the following Proposition [31 which 
provides a simpler version of Proposition [1] when the predator z does not exhibit density- 
dependence. 

Proposition 3 Suppose that the dynamics f is given by 

f(y,z,v,w) = I V y Vl ' ^ ] for all (y, z, v, w) G M 4 . (17) 
y zR z {y,w) J 

Assume that the function R y : M 3 — > M is decreasing in the control v and R z : M 2 — > M is 
increasing in the state variable y and continuously decreasing in the control variable w, and 
satisfies \im w ^ +QO R z (y, w) < 0. // 

—A > 1 and R z {y\ -r 
y> z 9 

the viability kernel associated with the dynamics f in ([17]) and the acceptable set D in (TTj 
is given by 

V(/, D) = | (y, z) | y > y\ z > z\ yR y (y, z, ^-) > y" \ . (19) 



R y (y\ z\ — ) > 1 and R z {y\ — ) > 1 , (U 



Proof. The proof follows the same lines as the one of Proposition [TJ The description of the 
viability kernel is simpler because of the following argument. The monotonicity assumptions on 
each component of function R z and condition R z (y 1 ', ^-) > 1 lead to 

y > y 9 and z > z° R z (y, —) > R z (y, > R z (y\ > 1 = Z - > - . 

z z 9 z 9 z 9 z 

Therefore, zR z (y, — ) > z^ for all y > y^ and z > z 9 . This explains why (|16p reduces to (|19p . □ 
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3.4 Viability kernel for a Lotka— Volterra system 

Consider the following discrete-time Lotka- Volterra system of equations with density-dependence 
in the prey 

y(t+l) = Ry(t)--y 2 (t)-ay(t)z(t)-v(t)y(t) , 

K (20) 

z(t + 1) = Lz{t) + (3y{t)z{t) - w{t)z{t) , 
where R > 1, < L < 1, a > 0, (3 > and k = ^^K, with K > the carrying 
capacity for prey. The dynamics / is as in ffTTl) . with R y (y, z,v) = R — — az — v and 
R z (y, w) = L + (3y — w . Proposition [3] gives the following Corollary 

Corollary 4 Consider the Lotka-Volterra predator-prey model ( 1201) . Whenever 



y b > 



l-L 



and z < 



R-l R(l-L) 



a 



any minimal catch thresholds Y b and Z b such that 



Y" <F b '* := y\R - % b - az b - 1 



R , 



Z' <Z b '* := z\L + j3y" - 1) 



(21a) 
(21b) 



satisfy (TT8I) and the viability kernel associated with the dynamics f in ( I20l) and the acceptable 
set D in ( 1141 is given by 



V(/, 



y>y\z b <z< 



a 



R 



k — y\ F b + y b 



(22) 



This Corollary has the following practical consequence. An initial state (y, z) such that 
y > y" and z > z b belongs to all the viability kernels V(/, D) associated to Z b < Z b '*, 
by (l2Tbl . and to Y b < min {y b >*, y(R - Ryjn - az) - y b }, by (l2Tij) and ([22), as long as this 
last quantity is nonnegative. In other words, if viably managed, the fishery could produce 
at least min {Y b '*,y(R — Ry/n — az) — y b } and Z b '*, while respecting biological thresholds 
y b and z b . We shall use this property in the following numerical application. 
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4 Numerical application to the hake— anchovy couple 
in the Peruvian upwelling ecosystem 



We provide a viability analysis of the hake-anchovy Peruvian fisheries between the years 
1971 and 1981. For this, we shall consider a discrete-time Lotka-Volterra model for the 
couple anchovy (prey y) and hake (predator z), for which the viability kernel has explicit 
description. The emphasis is not on developing a biological model, but rather on decision- 
making using such a model. 

The period between the years 1971 and 1981 has been chosen because the competition 
between the fishery and hake was reduced due to low anchovy catches, and because of the 
absence a strong warm event after El Nino 1972. We have 11 couples of biomasses, and 
the same for catches. The 5 parameters of the model are estimated minimizing a weighted 
residual squares sum function using a conjugate gradient method, with central derivatives. 
Estimated parameters and comparisons of observed and simulated biomasses are shown in 
Figure [H 

We consider values of y 9 = 7 000 000 t and z = 200 000 t for minimal biomass thresholds 



and values of Y b = 2 000 000 t and Z b = 5 000 t for minimal catch thresholds (IMARPE 



2000 



20041 ). Conditions ([TBI in Proposition [3] are satisfied with these values. Indeed, the 



expressions in (I21al) - fl21bj) give: 



Y b = 2 000 000 t < F b '* = 5 399 000 t 

(23) 

Z b = 5 000 t < Z b '* = 56 800 t . 

The viability kernel is depicted in Figure [2j The unique viable point (within the viability 
kernel) is the initial point. Thus, based upon this model, the fishery could have been managed 
- with appropriate viable controls - to produce catches above F b = 2 000 000 t and Z b = 
5 000 t, while ensuring biological conservation. What is more, due to the remark following 
Corollary H captures up to Y h = 5 399 000 t and Z b = 56 800 t were theoreticaly achievable 
in a sustainable way starting from year 1971. 
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Observed 

Lot k;i- Vol terra with 
ilmtiit.v-ik'jieiuli'neT in prey 




- Observed 
■ Lotka-Volterra with 
density- depends nee in prey 




(a) Anchovy (b) Hake 

Figure 1: Comparison of observed and simulated bio masses of anchovy and hake using 
a Lotka-Volterra model with density-dependence in the prey. Model parameters are R = 
2.25 year" 1 , L = 0.945 year" 1 , k = 67 113 xlO 3 t (K = 37 285 x 10 3 t), a = 1.220xlCT 6 t -1 , 
(3 = 4.845 x 1(T 8 t _1 . 

5 Conclusion 



Motivated by viable management of ecosystems, we have provided a general condition ensur- 
ing an explicit construction of viability kernels and have applied this to the viability analysis 
of generic ecosystem models with harvesting. 

Our results have then been applied to a Lotka-Volterra model using the anchovy-hake 
couple in the Peruvian upwelling ecosystem. We showed that, during the anchovy collapse, 
theoretically the fishery could have been viably managed to produce catches above the ex- 
pected threshold levels while ensuring biological conservation. 

It is interesting to notice that the kind of maximum sustainable yields V b '* = 5 399 000 t 
and Z h >* = 56 800 t provided by our approach in (1231) are comparable to the 4 250 000 t 
anchovy vie 



([PRODUCE 



and the 55 000 t hake yield, respectively, established for the year 2006 



2005 



2006b 



cj), or to the 5 000 00 t anchovy 



yield and the 35 000 t hake 



2006al 



2007a 



b). So, despite 



yield, respectively, established for the year 2007 ([PRODUCE] . 
simplicity 3 ] of the models considered, our approach may provide a mean of designing sustain- 



3 In addition to hake, there are other important predators of anchovy in the Peruvian upwelling ecosystem, 
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Figure 2: Viability kernel (in grey) for a Lotka-Volterra model with density-dependence in 
the prey in the predator-prey phase space (with y b = 7000 000 t, z° = 200 000 t, = 
2 000 000 t, = 5 000 t). The unique point within the viability kernel is the initial point, 
while all subsequent points in the trajectory are outside the state constraint set. 
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able yields from an ecosystem point of view. 

Thus, control and viability theory methods have allowed us to introduce ecosystem consid- 
erations, such as multispecies and multiobjectives, and have contributed to integrate the long 
term dynamics, which is generally not considered in c onventional fishery management. Notice 



that the World Summit on Sustainable Development (jGarcia. Zerbi. Aliaume. Chi, and Lasserre 



20031 ) encouraged the application of the ecosystem approach by 2010. 
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